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FOREWORD 


This paper was invited to provide a review of the portion of the registra- 
tion process relating to determining relative positions of reference and regis- 
trant Images at a set of locations within the images and embodiment of those 
positions in a reference-to-registrant mapping function. There was an unfor- 
tunately short time available between request and conference date, and as a 
result the review is not as thorough as the authors would have liked to have 
given, Ve are most familiar with tie processing of earth observations from 
Landsat, and have concentrated our review on the processing of Landsat data by 
the Master Data Processor (MDP) at the Goddard Space Flight Center (GSFC), the 
registration processor Installed at the Johnson Space Center (JSC), the regis- 
tration process used in the DAM (detection and mapping) Package in conjunction 
with a Cx>rps of Fjigineers project to map surface water, the registration process 
used by the LACIE processor ^t GSFC, and the sequential similarity detection 
algorithm (SSDA) used by IBM/ Houston in the evaluation of the LACIE processor 
imagery. Our review slights the body of information available from military 
researcli, in which shape recognition aud artificial intelligence play a promi- 
nent part. In that field, there is interest in finding objects that have 
changed position with respect to a fixed background, or monitoring changes in 
aspect for guidance, and there is often limited interest in spectral develop- 
ment- (Xir own field of concentration involves images in which considerable 
spectral development has taken place (observations are made throughout a crop 
year for agricultural applications), although for the mo«^t part the scene's 
basic geometrical sliape is unchanged. To some extent field boundaries change 
from year to year, but this is usually a small enough effect so as not to inter- 
fere with the registration process. In the presence of spectral development, 
one does not ordinarily correlate between images on the radiance measurement 
making up the parent image, but instead relies on a derivative feature (in our 
case, edges) that one hopes will be stable despite spectral development. We 
defer the discussion of extracted feature selection to other papers being pre- 
sented at this workshop, noting in passing that we are most familiar with the 
edge-detection algorithm implemented in both the and the JSC registration 
processor. This paper further concentrates on automated inter image correlation; 
the manual techniques involved in the DAM Package and other systems are a funda- 
mentally different technology. First, of course, the points of correspondence 
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between images are manually determined. Second, the interpretation Is normally 
done on radiance images transformed into visual presentations. And finally, 
the manual techniques are typically done on a point basis with only contextual 
involvement of the neighborhood around the point designated as being in 
correspondence. In the image comparisons discussed in this paper, measures 
of similarity are made over areal extents of subsets of the parent images, as 
opposed to points. 
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Introduction 


Interimage matching is the process of determining the geometric trans- 
formation required to conform spatially one image to another. In ^^rinciple, 
the parameters of that transformation are varied until some measure of "differ- 
ence** between the two images is minimized or some measure of "sameness** (e.g.* 
cross-correlation) is maximized. The number of such parameters to vary is 
fairly large (six for merely an affine transformation), and it is customary 
to attempt an a priori transformation reducing the complexity of the residual 
transformation or subdivide the image into small enough match zones (control 
points or patches) that a simple transformation (e.g., pure translation) is 
applicable, yet large enough to facilitate matching. In the latter case, a 
more complex mapping function is fit to the results (e.g., translation offsets) 
in all the patches. The methods reviewed have all chosen one or both of the 
above options, ranging from an a priori along-line correction for line-dependent 
effects (the **high-f requency correction'*) to a full sensor -to-geobase trans- 
formation with subsequent subdivision into a grid of match points. 

There is, of course, a correct geometric transformation to apply, but it 
is unknown (otherwise, an empirical image matching process would be unnecessary); 
thus, some sort of model must be assumed whose parameters can be solved for by 
correlation of offset-fitting. Commonly, a portion of the geometric model is 
established a priori based on external data such as preflight measurements. If 
that portion is incorporated in an a priori transformation, then the demand for 
fidelity in the overall model becomes an issue of tradeoffs between the a priori 
geometric correction and the residual correction to be determined by matching. 

For example, if one know the geometric transformation is affine with an unknown 
translation, one might let the Image matching function solve for the whole 
affine transformation, risking the introduction of errors into the affine co- 
efficients, or apply the affine transformation a priori and solve only for trans- 
lation. Thus, in addressing the various parts of image matching in the following 
sections, we must also consider their interaction with portions of the overall 
geometric correction process which otherwise! might be regarded as beyond the 
scope of this paper. The utility of an a priori transformation is also mani- 
fested in another way. Data gathered from prior registrations can be used to 
bias the a priori correction. Such **experience** data might be gatlieied by 
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analyzing ensembles of prior registrations such as associated acquisitions 
(different dates, same scene) or associated geometry (e,g,, different subscenes, 
same frame). Such experience data are then included in the data base. 

The following sections logically divide, according to the above precepts, 
into (1) a consideration of correlation techniques, (2) examination of other 
matching methods, (3) a discussion of determining the translational offset from 
the correlation, (4) a consideration of the residual geometric model and how to 
determine it, and (5) a summary of techniques going beyond the assumption implicit 
in the previous sections that match points have been established to allow correla- 
tion for translational offsets only. Throughout these discussion; we make refer- 
ence to the registration processor recently installed at JSC, note experimental 
results pertaining to that system as appropriate, and note general limitations 
and areas deserving further study. 

Image Correlation 

Suppose that an a priori correction has been applied or patches are defined 
small enough that any residual geometric error worthy of consideration is pure 
translation. Matching, then, consists of determining the translational offset of 
one subimage from another ''reference'* subimage corresponding to the same scene. 

The most common form of matching is some sort of cross-correlation technique. 

The "image" referiec to is the actual radiance image or a feature-space image 
such as an edge or gradient image. In principle the cross-correlation of two 
acquisitions of the same scene should resemble the autocorrelation displaced by 
the same amount as one image from the other. Since the peak of the autocorrela- 
tion function must occur at the origin, the peak of the idealized cross-correlation 
then measures the displacement. lu practice, the acquisitions differ because ol 
instrument, atmosphere, and other environmental noise and because the scene may 
have changed appearance due to seasonal, weather, or cultural changes. Thus, 
the cross-correlat ion only approximates the autocorrelation, and the peak may 
not be well-defined. Another argument for the cross-correlation peak measure 
is that it minimizes the sum-square-difference between the two images. Tliis 
measure of cross-correlation can be normalized (to be between -1 .uid 1) by two 
techniques summarized in Figure 1. The "template matching"^ alternative utilizes 
only the pixel values from the search area, whereas the "classical" alternative 
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utilizes both images. The pixel values Indicated are after subtraction of the 
image means. For binary edge images, whose pixel values are either 0 or 1, 
rigor fails in the following sense. Strictly speaking, the means of the O’s 
and l*s should be subtracted before summing; however, considerable gain in 
computation speed is achieved if the coincidences and numbers of I’s are 
summed directly. 

The correlation techniques just described are employed in several produc- 

2 3 

tion registration systems. The Master Data Processor (MDP) * , the GSFC baseline 

system for Landsats 2 and 3, basically performs the ’’classical** cross-correlation 

on the radiance image. As a f^^nctional equivalent, Fourier techniques are used; 

that is, the reference and search areas are transformed by the standard FFT 

algorithm, one of the transforms is multiplied by the complex conjugate of the 

other, and the result is inverse-transformed by FFT. A slight variance is 

effected in the denominator by subtracting the local mean (mean over portion of 

search area overlay ed by reference) from the search area's pixels, rather than 

the mean of the whole search area. The latter is used in the numerator in con- 

4 

formance with the classicai formula. Tlie GSFC LACIE Registration Processor 
employed the template matching algorithm on binary edge images, and the same 
scheme was used in an evaluation of that processor.^ Both classical and template 
matching algorithms are provided for binary edge correlation in the JSC Regis- 
tration Processor.^’ ^ The classical option has been chosen for production 
processing . 

The matching policy discussed above is based on minimizing the sura-square- 

difference. The sequential similarity detection algorithm (SSDA) is based on 

minimizing the sum-absolute-difference. However, rather than summing over all 

pixels, the SSDA selects at random a subset of pixels to sum. Summing stops 

when a present threshold is reached, and the number of pixels needed for that 

sum is noted. Then, the correlation maximum occurs at the same point as the 

maximum of the numbers of pixels required (because the sum-absolute-difference 

is smallest there requiring the most pixels to reach the threshold). The utility 

of the SSDA lies in its reduction in computation by utilizing only partial sums 

whose computational rigor is less, the smaller the correlation is. Thus, the 

bulk of the computation is devoted to correlation samples showing promise with 

9 

little effort wasted on low values. The SSDA was also used this time on the 
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radiance images, in the performance evaluation of the GSFC LACIE Registration 
Processor. The correlation was normalized by dividing the pixels of each image 
by their standard deviation before correlating^ 

The design of the SSDA illustretes a primary concern in standard cross«- 

correlation, viz diminishing the number of required computations. In general, 

methods to address that issue have focused on devoting the bulk of the computa- 

tions to the region around the appropriate extremum. Such a method, well-suited 

4 

to uinary edge correlation, was used in the GSFC LACIE Registration Processor 
and for coarse acquisition in the JSC Registration Processor. First, cross- 
correlation is done only for offsets every fourth row and every eighth column. 

The mean and standard deviation of the results are computed, and thrice the 
latter is added to the former to generate a rejection threshold. Next, an 
approximate correlation is computed at every offset but with only pixel values 
from every third row and third column included in the sums in Figure 1 (thereby 
reducing the computation by about a factor of nine). If, at a given displace- 
ment, the resulting approximate correlation value exceeds the rejection threshold, 
the sums are repeated using every pixel value. In binary edge correlation the 
peak is generally fairly sharp, and only a very small fraction of the correla- 
tion samples exceed the threshold, thereby sparing a considerable amount of 
computation. Inis method, used for a segment-level correlation (8x10 km por- 
tion of a Landsat MSS frame), proved very effective in the GSFC and JSC Processors. 

The foregoing discussion assumes the geometric difference between the com- 
pared patches is pure translation. This assumption may not be valid, which 
raises the question of how geometric distortion (other than translation) affects 
the cross-correlation function. The effect has been studied for standard cross- 
correlation (i.e., for minimizing the sum square 'ifference) where one patch is 
distorted linearly (rotation, scale change, shear distortion) from the other 
The results indicate that linear distortion effectively blurs the cross-correlation 
function by applying a running average over the ideal nondistorted counterpart, 
where the dimensions of the averaging filter are proportional to the degree of 
distortion. Tliis conceptual averaging does not apply to the noise present, so 
that the signal-to-noise ration (SNR) is effectively reduced. Also, the area 
around the peak is effectively flattened somewhat, making the peak search less 
accurate, and the peak-to-background (or peak-to-sidelobe) ratio is distorted 
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(these subjects are discussed In a subsequent section). The magnitude of these 
effects is greater, the greater the patch size, because the geometric distortion 
becomes larger. Thus, choosing a smaller patch decreases the effects of geo- 
metric distortion, but the SNR is also decreased with fewer pixels. Some 
compromise is generally required, as is discussed In a later section. 

The effects of geometric distortion on binary edge correlation could be 
drastic in certain cases. For example, a scale difference might prevent over- 
laying edge features on opposite sides of a patch at the same time. Thus, rather 
than smearing the correlation function, as described above, it could divide into 
several peaks. This kind of problem has not been identified in typical patch 
correlation, but it has been suspected of degrading the full sample-segient 
correlation used in the GSFC LACIE and JSC Processors. One solution to the 
problem (at least for coarse correlation) might be to resample the images to a 
coai^er grid by use of a low-pass or median filter before edge detection. Then, 
the geometric distortions are reduced relative to the pixel spacing or, equiva- 
lently, effective edge "thickness . 

Other Techniques for Image Offset Matching 

Sever< 1 other matching techniques, though not yet implemented in large- 

scale registration, deserve mention as potentially applicable. The Cluster 

13 

Reward Algorithm (CRA) con itually analyzes the bivariate histogram of the 
two images at each displacement. A measure is established characterizing the 
definition of pattern in the histogram, e.g., a measure of clustering. At the 
offset denoting an image match the histogram should show a relatively high 
image-to-image correlation by exhibiting clusters. The technique has been 
applied to several sample-segnent-size (about 8x10 km) scenes with notable 
success. 

14 15 

The point-matching technique * , though designed primarily for target 

arrays such as star fields, might find application to binary edge images. The 
method minimizes a geometric distance measure between points (say, edge pixels 
in the two images) rather than ninlmizlng radiometric or feature differences. 

A variation^^ of the method ascribes weights to the point-pairs according to 
how well they associate in the matching. This technique is iterative with first 
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a matching, then a comparison of the distance spanned by each point-pair to the 
estimated common displacement for assigning a weight, then a weighted match, etc* 

The use of normalized cross-correlation for minimization of um-square- 
differences or sum-absolute-differences described in the previous section can 
also be viewed as an application of maximum-likelihood classification in which 
the measurement sets are distinguished by translational offset (analogous to 
pixel number in spectral classification), the channels of information are the 
pixels in the search overlay area (analogous to the bands in spectral classifi- 
cation), and the classes are Match and Non-match (analogous to classes, e.g., 
crop species, in spectral classification). In the spatial case, however, it is 
known a priori that only one measurement set should be classified as a Match (if 
the match scene is unique) corresponding to the single correct translational off- 
set. Analogous to spectral classif icatic a Bayes technique has been applied^^ 
to spatial matching. The probability of misacquisition (match lock-on to the 
wrong point), expressed as the Bayes risk is minimized by maximizing the a 
posteriori probability of a correct match (i.e., the probability that there is 
a match at an offset, given the particular pixel values in the overlay region 
for that offset). The a posteriori probability is estimated in principle by 
least squares by assuming it is a linear function of the pixel values. In fact, 
that probability should exceed 1/2 only at the match offset which provides a 
direct means for identifying a probable misacquisition (a posteriori probability 
never exceeds 1/2). The method is conveniently implemented by Fourier trans- 
formation techniques. The a posteriori probability, though assumed linear in 
the pixel values above, can be expressed as a linear combination of any functions 
of the measurements, as deemed appropriate to the application. The case of a 
second-order relation has been investigated^^, and trials on synthesized data 
indicate improvement over the maximum likelihood alternative offered by cross- 
correlation. 

The methods described in this section should be considered as new candidates 
for interimage matching. They must be integrated with the various feature selec- 
tion methods (e.g., the Cluster Reward Algorithm makes no sense with binary edge 
images, while the point pattern matching technique is especially suited to such 
images) . The methods should be tried on a representative set of sensor data and 
evaluated for applicability to different sensor types, scene types, season, 
sensor-to-sensor differences, etc. 



Offset Determination 


The cross-correlation tecf*iiiques for image matching presun»es that the posi- 
tion of the maximum ot minimum cf the correlations or difference functic.. charac- 
terizes the translational offset between cbe two images. The true offset might 
be in between pixel centers so that the correlation peak lies between correlation 
sauries. To achieve subpixel accuracy in a given patch, some form of interpola- 
tion is needed. Several alternatives have been dc /eloped, as we discuss in 
this section. These techniques lend tberaselves to estimating the accuracy of 
peak location and to evaluating certain measures for establishing pass/fail 
status, and vre address this subject also. 

As an alternative to peak inte eolation for subpixel accuracy, a number of 
patches can be definid, the peak can be determined to the nearest sample, and 
the subsequent mapping function fit to the numerous peak offsets can be relied 
upon to yield subpixel accuracy. This approach is predicated on the assumption 
that a random position error is introduced by choosing a sample position rather 
than interpolating. Tliis method was used in matching images produced by the 

GSFC LACIE Processor for the purpose of evaluating the performance cf that 
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processor*^* . Second«»ry peaks were also identified and compared t^ the primary 
peak. A measure of decline away from the peak was established in terms of 
averages of correlation samples in concentric rings extending out from the peak. 
If the peak was not strong enough relative to the secondaries or the decline was 
too shallow, that zone was rejected in the mapping function fit. 

A straightforward means of interpolating for the peak consists of fitting 
a surface function to the cross-correlation samples in the vicinity of the -^eak 

and evaluating the peak location analytically. The KDP and JSC Registration 
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Processor * use a bivariate polynomial (i.e., containing terms x y , MN ^ poly- 
nomial order). The MDP uses a fourth-order polynomial on a 5x5 neighborhood 
around the peak and also evaluates the curvature at the peak as a measure of 
the breadth of correlation peak. The smaller the curvature is, the larger 
is the uncertainty in locating the peak. The MDP uses the minimum curvature 

on the surface and the height of the peak as rejection criteria. 
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The JSC Registration Processor allows any order polynomial up to and in- 
cluding fourth order, but the latter has been used mainly, with the lower orders 
18 

being studied for applicability* The fit neighborhood size can also be 
18 

varied , lut a 5x5 has been used in prod»iction* The peak is located by the 
use of the MDP algorithm, and the curvature matrix is used, in lieu of the 
curvature itself, in a variatimal approach to express the peak uncertainty 
covariance in terms of the covariance of uncertainty in the polynomial fit co- 
erificients. Since the bivariate polynomial is fit by least squares, the coeffi- 
cient uncertainty covariance can be estimated from the fit residuals. Often, 
relatively few points are fit, compared to the number of coefficients solved 
for (25 points for a 5x5 neighborhood, as compared to 15 coefficients for a 
fourth-order pol 3 Tiomial) , and little confidence can be placed in the uncertainty 
estimate. To circumvent this problem, the user may supply an a priori fit vari- 
ance of the fit uncertainty rather than determining it from the residuals. The 
a priori variance can be established by tests on representative data, whereby 

the a priori variance can be adjusted to make the estimated peak uncertainty 

18 

agree with the observed peak displacement errors. Preliminary tests estab- 
lished a value for early production in the JSC Registration Processor. A peak- 
to-background ratio is also established by subtracting from the correlation peak 
value the mean of samples away from the peak and then dividing by the standard 
deviation of those samples. The subtraction was included to compensate for the 
fact that the edge image means are not subtracted out in binary edge correlation. 
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Rather than establishing the orvler oi' tho bivariate polynomial a priori, 
saccoasive orders could be tried and the RMS residuals compared, Vrtien that 
error falls significantly and then starts tapering off, whatever order has 
been used at that point is adopted as the fit. As an attempt to sidestep tho 
Issue that a g^x>d model of the oross^orrelat lou is tu^t really at hand, this 
approach follows tho rationale that a good reprosentat ion lies somewhere be^ 
tween a simple model and a fit, with iu> remaining orr^^r, to every point, as a 
sufficiently high^ordor polynomial would do. Choosing an intermediate order, 
as outlined, attempts to balance a lack of understanding of tlie general form 
of the cross-correlat ion function with the knowledge that errors do exist In 
tho data which makes forcing a perfect lit unwairantod. Tliis pixd^lem of model 
uncovt dinties is eucouatored again in tho section dealing with the mapping 
fuactioit. 

FunctivM\s other than a bivariate polynomial can be fit to tho cross- 

correlation samples. A bivariate gausstan function Is an obvious choice; how* 

l'^ 

ever, its use has not come to our attention. One registration application 
chose an elliptical cone. Tlio orientation of the cone was also a solution 
variable. Although surface fitting !\as boon tlio Cv>i’rnonost method of peak 
Interpolat Ion, it is bv no means universally accepted. A principal criticism 
with that method involves the sens it Ivltv of the inters;imple peak local ivm to 
the particular form of th*e fitting function. 

As an alternative to surface fittit\g, the peak Cxiu be assumed to lie at 

'0 

the centroid of the correlation neighborhood ' . T!\at is, tl\c peak location 

is computed as the weighted sum ol neighboring sample locatloixs, with tho 

correlation v; Uies serving as weights. Analogous to the peak uncertainty or 

surtaco curvature meat ioned earlier, a corresponding measure of the breadth 

:i 

ot the peak can be computed as the moment ot gyrativm or simply as the 
second moment of tho correlat ion distrihutloa ijust as tho centroid Is t!ie 
first moment). Unfortunately , Htllo more can be said about this metliod at 
this t imo . 1 1 wou I d be i n to res ting to c ompa i e t he c on t ro 1 d auvl su v lace t i t 

me t liods tor a r ep resen t a t i v e set o f cv' r re I a t i on samp I e a r i y s . 



Another method of peak location takes advantage of cross-correlation by 
Fourier techniques. The offset in the spatial domain corr^^sponds to a non- 
vanishing phase in the frequency domain. In fact, if the cross-correlation 
function were symmetric about its peak, the phase function would directly 
specify the offset (since phase introduced by non-symmetry would vanish). 

'70 

One approach*''^ transforms the phase portion of the correlation transform back 
to the spatial domain wherein the peak is located. The spatial result was 
noted to resemble closely a delta function so that subpixel peak location 
seems viable. The peak location was facilitated by using the inverse transform 
operation directly to compute a finer grid of spatial samples around the indi- 
cated peak than would normally be obtained from Fast Fourier Transform (FFT) 
algorithms. If indeed the transformed phase function is consistently nearly 
a translated delta function, then in the frequency domain the phase should be 
proportional to the spatial displacement. Thus, the phase vs frequency points 
can be fit by a straight line forced to pass through the origin, and the slope 

gives the translational displacement. However, somehow the possibility of non- 

18 

symmetry may need to be accounted for. Preliminary results have shown that 
binary edge correlations, for example, do exhibit nons\mimetry. 

This section closes by returning to the issue of estimating peak location 
uncertainty or sharpness. The schemes described above can be categorized 
basically into curvature at peak, second moment of correlation samples around 
peak, and rate of decline away from peak. Any of these sharpness measures can 
be used to compare with a rejection threshold or to construct a weight for use 
in the mapping function fit. In the latter case the weight would be lower, the 
broader the peak lobe is. The weight is given directly by a peak uncertainty 
estimate (as the inverse of the uncertainty covariance). As mentioned earlier, 
a fixed a priori surface fit variance should probably alwavs be used to estimate 
the peak uncertainty when fitting bivariate polynomials so that effectively its 
, variance is directly related to the curvature. The sharpness measures can 
also be computed for the autocorrelation function to give an idea of the intrin- 
sic clarity of the image. This application has been implemented as well in the 
JSC Registration Processor. 
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Residual Mapping Model 


If the approach of the previous sections has been adopted, viz. local 
translational mismatches are determined at an array of control points, some 
means is needed to distribute the corrections for those mismatches over the 
regions in between the control points. Toward that end, a geometric transformation 
function is to be lormulated in terms of the measured displacements . The .uiture 
of that function depends upon one’s confidence in the measured mismatches versus 
one’s understanding of the geometrical or physical processes causing the mismatch. 
If the control point translations are expected to contain errors, and the nature 
of the geometric distortions is well understood, then some sort of error mini- 
raization or maximum likelihood estimation technique is utilized. On the other 
hand, if the data are highly trusted and the model is unknown, some form of 
’’rubber sheeting" technique is employed (analogous to curve-smoothing in one 
dimension). These alteri'iatives are addressed below. Tl^e case that the model 
is understood is discussed first, with a breakdown into a geometry-based model 
(e.g., platform attitude error, trajectory error, etc.) and a sensor-based 
model (e.g., scan irregularities, band-to-band offsets, etc.). Then, the rubber 
sheet approach is considered. If uncertainty is associated with both the measure- 
ments and the model the attack is not so clear. This situation is considered 
briefly along with the situation t!iat some measure of cotifidence in the model 
is available and the control point data are being used to improve it. It is 
well established that cross-correlat ion will occasionallv result in an erroneous 
translational offset, and if such an error goes undetected it can ruin the mapping 
fit. The section ends with a description of several methods for identifying 
and handling such outliers. 

Geometric errors are generally modelod as time varying satellite attitude 

3 

and altitude errors. The MDP assumes tor MSS processing that the I^ndsat yaw, 
pitch, and roll can be modeled as third-order functions of time over the dimen- 
sions of a double frame (about 340 km of downrange) . Similarly, the altitude 
variations are modeled as linear in tune. Tlie resulting 14 coeff jcients are 
solved for by weighted least squares utilizing botli control point offsets and 
the loss-precise attitudc/alt itude data available from Uindsat *s attitude 
measurement system (AMS) and ground tracking. T\\c weights are set up a priori 
based on past experience in control point location accuracies and in AMS and 
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ground tracking accuracies Hence, the control point weights are considerably 
larger than ^ae others so that if many control point offsets are utilized the 
solution is essentially based on them* However, if few or none are available 
the AMS and ground tracking data do conttJuDute appropriately to the solution. 

To handle the fact that the AMS data are poor for estimating the constant por~ 
tion of the model (biases), a scheme is implemented to allow devoting the con- 
trol point data to estimating the attitude biases, while the other measurements 
are utilized for rates and higher-order coefficients. Nominally, this scheme 
is employed only if a very few control points are available. Also, the attitude 
biases from a good solution (one with control points) are carried over to the 
next frame i‘ that frame has insufficient control points (nominally, less than 
two). A weight associated with that carry-over is propagated in a manner that 
makes it decay with time, in order to model the growth in uncertainty of the 
biases as the satellite moves further away from the estimation point. This 
whole approach is predicated on the measurement errors being nornu.lly-distributed, 
an assumption which breaks down if one or several of the control point offsets 
are erroneous. To overcome that problem, the MDP has been tested with the 

threshold number of control points for triggering the bias estimation scheme 

23 

raised to 15 from the nominal value of unity. Hopefully, such a large number 
of control points will contain enough good offsets to overwhelm any erroneous 
ones. Another approach to eliminating the outliers is discussed later. 

24 

A TRW study has adopted the same model as the MDP except for assuming a 
constant altitude deviation. The coefficients were evaluated by Kalman filter- 
ing. ERIM’s geometric correction process also models the geometric errors as 
yaw, pitch, and roll, with pitch and roll assumed linear and yaw assumed con- 
stant over the area of interest. 

The GSFC LACIE processor models the geometric error (after the more-compli- 
cated a priori transformation) as pure translation ever the dimensions of its 
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8 X 10-km sample segments. Performance evaluation * of that system indicated 
thac, although specifications were met, there was a residual error that appeared 
to be affine. Thus, the JSC Registration Processor has adopted an affine model, 

A weighted Icjast squares solution is made in which the weights are the reciprocals 
of the peak uncertainty standard deviations estimated from the correlation surface 
fit. At user option, default weights can be substituted if the correllation 
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surface fit fails in its estimation process, or weights can be eliminated alto- 

26 

gether* A preliminary evaluation of the JSC processor has indicated that an 
affine model is adequate; that is, observed residual distortion does not show 
any systematic effect associated with an incomplete geometrical description. 

The success of the affine model lies in the fact that the JSC Registrar ioi* 
Processor is only removii g residual distortion left by the MDP. The a priori 
geometric correction has taken care of map projections, sensor geometry, and 
sensor line effects as characterized by the MDP data (a priori modeling con- 
stants and attitude/altituJe solution). 

Sensor effects can also be modeled and solved for by use of the acquired 
scenery. Alternatively, they can be assumed static, determined from a priori 
data, such as laboratory or flight tests, and then used in the a priori geometric 
correction. Landsat MSS effects include scan nonuniformity, band-to-banu offset, 
**staircasiug” caused by simultaneity of six scan lines, and line-by-lin«> offsets 
due to sampling delays and changes in scan speed. The MDP, GSFC LACIE Processor, 
and JSC Registration Processor all assume a static model for those effects. All 
effects are assumed constant, except for the scan speed and scan nonuniformity. 

The latter is modeled as a third^ordcr polynomial in position with the line (i.e., 
sample number). The coefficients were determined prior to implementation in the 
registration systems. The scan speed is assumed to have a fixed relation to the 
number of samples actually obtained in a scan line (line length majority) - Al- 
though this model seems fairly good in normal circumstances, it breaks douTi for 
Landsat 3 when the MSS sample initiation fails (the "late line start") because 
the correct line length majority Is not available. A model has been installed 
in the GSFC processing line to estimate the line length majority from the number 
of samples in the partial line. The resulting line length majority is used by 
both the MDP and the JSC Registration Processor. An inaccurate line length 
majority is manifested as a line jitter over portions of the image, and iTieed 
a jitter is sometimes apparent as irregular field boundaries in some Landsat 3 
agricultural scenes (especially at lower lattitudes where field boundaries tend 
to be parallel to lines and columns of pixels). An error of one in the line 
length majority will result in a misregistration of up to one pixel or about 58 m 
in MSS imagery. The possibility of variations in scan nonuniformity can also not 
be ruled out as the cause. These anomalies merit further investigation. Tech- 
niques for line-to-line registration to straighten field boundaries should be 
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considered, if for no other purpose than clarifying operational scan 
characteristics . 

The sensor geometrical model for the Laudsat RBV is considerably simplified 
by the presence of reseau marks on the face of the Vidicon tube* Also, the 
satellite attitude model is simplified by the fact that a whole frame is 
gathered at one instance, hence, for one value of roll, pitch, and yaw. 

The rubber sheet approach to geometrical oodeling assumes the form of the 
model is poorly understood. Whatever is the correct model, it should be possi- 
ble to express it as a bivariate polynomial (e.g., assume a Taylor series). 

Thus, a common rubber sheet technique fits a polynomial to the array of control 

point offsets. The JSC Registration Processor and the GSFC Digital Image Recti- 
27 

fication System can utilize polynomials, the former to fourth order and the 
latter to second order. The principal concern is generally what order to use, 
although omission of certain terms may also be considered (such as sample- 
dependent terms if distortions seem to be from line to line). A sufficiently 
high order will fit all control points with dubious resax ts in between. If the 
points are known rigorously to contain no error, then fitting all the points is 
reasonable, and the problem becomes one of selecting an interpolating function 
with suitable behavior between the points. Gener^^l , the control points are 
assumed to have errors so that fitting them rigorously at the risk of inter- 
point error is not particularly suitable. A compromise is often achieved by 
trying a progression of successively higher-order polynomials and tracking the 
decrease of the error, expressed as the RMS of fit residuals. Assuming the 
idea that lower order is better for inter-point behavior, the fit is chosen for 
which the error has dropped substantially from lowest order but for which little 
error reduction is apparent at higher orders. 

Another rubber sheet approach effectively spatial-f il'^ers the ’’sampled** 
offsets at the control points, thereby distributing the offsets over areas 
between the control points. This approach parallels the reconstruction of a 
two-dimensional analog signal from its samples by cubic or sine convolution. 

The filter might be derived from maximum likelihood considerations such as a 
Kalman filter. This approach might be reduced to any desirable scale by reducing 
the size of control points while increasing their density so that a fairly 
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high-frequency matching process could be designed. When sensor resolution be- 
comes good enough, atmospheric turbulence will be discernible as high-frequency 
distortion for which sucli a high-frequency geometry-matching filter would be 
appropriate . 

Sometimes, the model may be partially known before control point informa- 
tion is acquired; that is, a priori values for the model parameters are available 
with an associated covariance characterizing the level of confidence. Then, the 
control point offsets provide additional information which can be used to improve 
the estimates of the model parameter. A Bayes approach is utilized. Simply 
stated, a least squares solution is performed on the control point offsets, and 
a weighted average of the results with the a priori estimate is computed. The 
weights are the inverses of the respective covariances normalized by matrix pre- 
multiplication by the inverse of the sum of the covariance inverses. This 
approach has been implemented in the JSC Registration Processor to account for 
the fact that it follows the MDP. Since the MDP has already geometrically 
corrected the data, the a priori values for the residual mapping function are 
zero. The covariance of the MDP*s attitude/altitude solution is provided, and, 
after transformation to the covariance of a local (8 x lO-km sample segment) 
affine geometrical model, it serves as the a priori covariance. The JSC Proc- 
essor performs its own least squares solution and utilizes its and the a priori 
covariances to average its solution with zero (the MDP*s value). The appropri- 
ate expressions are summarized in Figure 2. Effectively, the result down-weights 
the JSC Processor's solution according to the size of the MDP covariance, with 
greater •attenuation the smaller the MDP covariance. It must be borne in mind 
that all this theory is predicated on the MDP's and JSC Processor's error sources 
being normally-distributed. A crude compensation (shrinkage) is provided, as 
shown in Figure 2, in case they are not normally-distributed. 

Unfoi Innately, control point correlations tend to either work reasonably 

well, or they are very bad. Thus, unless false fix detection is very good, 

erroneous control point offsets will be interspersed among the good ones. One 

straightforward approach to eliminating the outliers simply comp, res the residuals 

after the least squares fit to a prese*^ thrc:shoid, or a threshold scaled by the 

estimated covariance. Thio technique was used in the GSFC LACIE Processor eval- 
5 9 

nation • . Control points with residuals greater thin three standard deviations 
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were discarded, and the least squares solution was repeated. The new residuals 
were again tested, and the procedure was repeated until no more failures were 
noted « 

Another approach performs a least squares solution with one control point 
CMnitted. Similarly, another solution is made with the first control point re- 
turned but a second omitted, and the process is repeated to yield as many solu- 
tions as there are control points. The post-solution residual for each omitted 
control point may then be examined for rejection, or a weighted combination of 
all the solutions may be made, with the weights related to the post-solution 
residuals. 

28 

Finally, a robust estimation technique can be utilized in which the con- 
trol point error is assumed to come, not from a normal distribution, but from a 
distribution with long tails which characterize the possibility of ’’outliers.** 
Just as least squares is a maximum likelihood estimator for a normal distribu- 
tion, maximum likelihood estimators can be formed for other longer-tailed 

28 

distributions. A number of such solution techniques exist. They down-weight 
measurements with large errors, while sacrificing a small amount of efficiency 
(i.e., not being quite as good as least squares) when the measurement errors 
really are normally-distributed. To our knowledge, none of these techniques 
have been explored for image matching. 

Sizing and Placement of Control Points 

'fhe method of inter-image comparisons reviewed here involve the distributing 
of ^ number of patches within each image (the reference image and the registrant 
image). The tches are chosen small enough that the inter-image shifts within 
them raiv be considered to be purely translational; that is, scale, skew, and 
other interlnage distortions have negligible effect when compared to the trans- 
lation, viewed on the scale of the patch. Translation we will regard as zeroth 
order, amounting to a constant bias in the functions that relate coordinate 
values in one image to coordinate values in the other. Similarly, scale, rota- 
tion, and skew are first order distortions, and keystoning, etc., are higher 
order. In the zeroth order, coordinates reference image and in the 

registrant image are related 
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(Xl-Yi) = , 


while to first order they are related 

(Xi-Yi) = +b^ . 

So, it is seen that the above description of the order corresponds to the order 
of the polynomial in the equations relating the coordinates. The essence of 
registration is to specify the mapping between the reference and registrant 
images, which is done in either "rubber-sheeting” form or by modeling. Ttie 
latter uses rigid descriptions of the sensor geometry to permit an economy in 
describing the interimage coordinate relationship. Sc, rather than describing 
the coordinate relationship in the form above, coefficients are often drawn 
that relate the coordinates through platform ephemerides and known sensor 
geometries. This section interrelates with the previous one via economic and 
accuracy considerations; an optimum mapping method will require the fewest con- 
trol patches to establish the desired degree of accuracy, since there is a high 
cost per correlation patch. 

Number, location, and size of the patches used for control points is a 
very important consideration in the design of a registration system, from both 
economic and accuracy standpoints. A large part of the registration computation 
load is proportional to the number of the patches and to the second power of the 
patch dimension, so it behooves one to minimize the number and individual size 
of the patches. At the same time the interaction of the spatial distribution of 
the patches with the remapping model will affect the accuracy of the registra- 
tion; so the location of the population of patches is a consideration to be com- 
bined with their number. One patch gives a single estimation of the translational 
correspondence between the images in the vicinity of the patch. If it be assumed 
that errors are not systematic in the estimation of that translational corre- 
spondence, more patches (in a neighborhood) are better in that the errors will 
tend to cancel. However, although there is no conceptual difficulty with having 
patches overlap, there is a limit to the number of functionally independent 
; itches that can be drawn from any image. That number puts an upper limit to 
the tendency for growth to larger number of patches in order to have statistical 
cancellation of errors. 
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There is a tradeoff in the determination of the optimuai size of a correla-- 
tion patch that is dependent on the accuracy with which the a priori corrections 
can account for the encountered misregistration. Consider first a purely trans- 
lational misregistration. As the size of the patch gets larger, there is more 
picture structure within it; although the value of the correlation peak is not 
affected, at translational locations off the peak there will be (generally) 
lower cross-correlation values due to the larger amount of structure within the 
window. The result is an enhancement of the peak-to-background ratio (the value 
of the peak correlation divided by the correlation in the general vicinity), 
enabling increased sensitivity in detecting a peak and perhaps increased accura- 
cy in determining the location of the peak. 

Itow, however, consider the ramifications of first-order departures from the 
condition in which there are only translational offsets between the images being 
registered (i.e., suppose there are skew and scale differences as well as trans- 
lation). In the position in which the true centers of the patches being corre- 
lated are coincident (that is, the position one would want to find as a result 
of the cross-correlation process), picture structure towards the edges of the 
patches will be relatively displaced by scale and skew differences. In contrast 
to the translation-only situation, the correlation in this case is reduced as 
the size of the image is increased. But the presence of effects other than 
translational offsets does not mandate the use of the smallest possible patch — 
one pixel! — because the correlation is at its noisiest there. For a binary ex- 
tracted feature on which correlation is being done, the correlation takes on 
only the values zero or one for the one-pixel patch, and it would be impossible 
to locate a peak. So one determines an optimum patch size by considering a 
reasonable envelope for the distortions (beyond translational) that one expects, 
and works out the spatial range over which they would cause significant displace- 
ments for the translationally-correct position. 

Where appropriate by necessity of increased registration accuracy, a re- 
entrant technique is possible. In that technique, patches of the original 
images are extracted and run through the correlation process. The resulting 
remapping for the whole image is then used to reextract the correlation the 
correlation patches, which can now be at least first-order corrected for the 
modeled distortions. Tlie patches ar<^ extracted from tue input imagery by 
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resampling onto a grid derived fiom the first registration process, and conse- 
quently larger patches, with potential for increased accuracy, may be drawn. 
The cost of such a reentrant technique is not necessarily to double the compu- 
tation time, inasmuch as a coarser registration job is enabled for the first 
pass. 


As mentioned earlier, the strategy for location of the correlation patches 
is important. In the case of the MDP, there is a scattering of patches amount- 
ing to approximately 0.1% of the scanned area (say 10 patches of 32 x 32 pixels 
within a frame of 3596x2983 pixels). In order for this small a correlated 
area to control the overall registration accuracy to an acceptable level, the 
patches must be lo. ited efficiently so as to exercise reasonable control over 
the entire image. Further, the modeling of the scanner behavior betweei control 
points becomes of paramount importance. In the other extreme, the JSC Registra- 
tion Processor actually uses something over 100% of the area of the registrant 
image, when it is taken into account that the patches overlap slightly (about 
7 pixels/line). 

The registration of two images logically breaks into two philosophically 
different kinds of operations — the determination of the mapping of the coordi- 
nate grid of the reference image into the coordinates of the registrant image 
onto that remapped grid of the reference image. The mapping process is guided 
by the selection of patches within which interimage comparisons give the local 
estimate of the mapping (typically just the translational characteristics are 
determined on an area small enough that translation is the dominant effect over 
the size of the patch). Since considerable computation effort goes into the 
interimage comparisons done on the patches, and since the accuracy of the over- 
all registration defends linearly on the accuracy with which a single patch is 
compared between images, and also since the model, usually driven by the trans- 
lational offsets determined at the patches, will be more tightly and accurately 
determined with a "sufficient** number and distribution uf patches, it is very 
important that the following qualities are achieved in the location of the 
control patches. 

1. The comparison method allows an accuracy of interimage comparison at least 
as good as the registration accuracy desired for the whole process. Other- 
wise one must use a superfluity of control patches and hope tnat they 
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contribute Independent measurements^ with error cancellation, to the re- 
mapping model. 

2. The model being driven by the control patches Is sufflci^^ : ...y highly de- 
tailed to account for the entire behavior of the relative geometry In 
between control points. 

3. Th3 locations of the control points are chosen to Interact sensitively 
with the modeling geometry. For example, if a sinusoidal component to the 
interimage geometry is hypothesized, one would not want control pc-^nts (the 
centers of the patches) placed at uniform one-period separations within the 
image. If they were, aliasing would have the result of a sinusoidal rela- 
tive gecxaetry look the same as a translational offset. 

This section discusses some aspects of location of control points (taken as the 
centers of their control patches). 

Let the registration of Image 1 to Image 2 be defined as the mappxng of 
the coordinates (i,j) in Image 1 to tne coordinates (k,"') in Image 2 (i. • , we 
will do only the geometric portion of the problem) . 


i k -> 


j 

IMAGE 

1 

IMAGE 


1 

4- 

2 


Suppose the "true relationship is given by and F, 


i = Fj^Ck.l) 


j = F^CkA); 


let them be estimated as a result of the intertmage correlation and modeling by 


fl = 


f =* F c 
2 2 


Then the root-mean-square registration, error (RMSE) is 


RMSE 


■\/ 15 r <F,-£2)2ldA 

'' image 
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where Is the Integral over the area of the Image, and various assumptions 

are included regarding orthogonality of the coordinate axes, samene s of scale 
on those axes, itc . Further, suppose the mapping f^ and f^ to be derived accord- 
ing to some sensor/platform model, or ac< ording to other limitations on Che com- 
plexity such as a mapping of no higher order than affine, and that the mapping 
and f^ is determined by a finite number of control pcnts within the Images. 

If the images are not Initially pathologically mlsregistered, we may reasonably 
express the mapping by beginning with the simplest of mappings, 

f^(k,l) - k , f2(k,l) * 1 

and then include the deviations g^, 

fj^Ck,!) « k + g^(k,l), f2 = 1 + g2(k,l). 

This slightly complicates the express on for RMSE but facilitat.s the analysis 
of sensitivity with respect to control point location. Translational offsets 
are wr-'.tten 

gja.l) = 

gjtt.l) - 

affine mappings have the general form 

g^(k,l) - + ^2^ 

g.(k,l) - b + b.k + b.l ; 

2 0x2 

keystoning is modelable by the inclusion of a cross term, 

g^(k,l) = ^3^ 

g^Ck,!) = b^ + b^k + b^l b^kl , 

etc • 

We consider a control point established between the images. Let ^(k^l) 
be the derived estimate of the *true** rffset between the images, and let I 
be the estimated covariance of the estimate t. 
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By a method unspecified here (so as co preserve generality), the set of control 
points {km, Im, ^m, Em} is then used to generate the mapp^^g tunctious and g^* 
Let us make some simplifying assumptions: 

1. The method of obtaining g^ and g^ will work with any number of control points 
(e.g., no control points gives g^ and identically zero, one control point 
allows a translational solution, etc*). 


2* The method of obtaining g^ and g^ is tractable to linearization. 

3. A priori estimates of I are available (they might be in the form 

^ap ■ "ap ‘ ‘ 

4. The a priori estimates of E are independent of location. 

5. The a priori estimates of g^ and g^ are zero. (Otherwise the image could 
be adjusted according to the a priori estimates so that this would be true.) 


We are now in position to find the sensitivity of the mapping functions and g^ 
to the location of a single control point, that sensitivity being a function of 
location of control po nt in the reference image and of location in the registrant 
image. The sensitivity can be converted into an rms value by integration over 
the registrant image, as will be seen, and thus an ideal location of the first 
control point is specified. 


i k 

Defining y_ = (j), x ” ^ have used the set {x, E} to produce the 

(vector) fu* . Ixon such that 

= g^(x) -h k. 

The fact that the variation in, say, the first coordinate of v is given by 


3gi 

= IT" ^ 


6x^ 


leads us to the vectoi variation in 
6^ = ^ ox , 

in which J is the Jacobian matrix 
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ORIGINAL PAGE IS 
OF POOR QUALITY 



The covariance In is calculated from the covariance in x, as follows. 

* ^ ^ si ^ * 

T 

Taking expectation values, and assuming that the expectation of ^ 6^ is as 
estimated from the considerations of the cross-correlation surface for the point, 

» E{((5y) (6y)^} « £{(6x) (5x)^) ^ 

S “ ] I J 

y ^ X ~ 

which will be a non-constant function of ^ inasmuch as the partials forming 
are not constant. The trace of gives, as a function of v, the squared regis- 
tration uncertainty. We can average over the image and take the square root. 

t 

RMSE = V / ff tr E dA 

V area . y 

^ image 

Alternatively the trace can be taken after integration (the operations conmunite) 

r — 

RMSE - \ / tr ( . // E dA) . 

y area image y 

Note that the RMSE so calculated is a function of x* the location of the first 
control point. If a priori values for and ^ (x) are used (the latter being 
zero), we have an a priori estimate of the sensitivity, measured in terms of 
coordinates in the image, as a function of location of the first control point. 
Fixing the location of the first control point at the location that minimizes 
RMSE, then an exactly similar procedure leads to the location of the best place 
to position the second control point. The RMSE should mono tonic ally decline as 
the number of optimally-placed control points increases, and the system designer 
can quit adding control points when the anticipated accuracy has gotten to 
requirement level. 
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It is true that the placement of n control points by this technique will 
(probably) not produce as good accuracy as a points placed optimally as an en- 
semble. The algorithm for determining that optimal ensemble position, however, 
is a substantially more complicated one. The following relaxation technique is 
such an algorithm. The n points are placed according to some scheme in the image, 
and RMSE is calculated in the general manner as described above. The position 
of the first point is permitted to vary about its original location, and RMSE is 
observed; the location is altered until a minimum of RMSE against first-point- 
location is found. Then, the second is allowed to vary, and its minimum location 
found. Similarly through the set of n control points, and the procedure is 
iterated for all the control points again until each point is at a minimum. But, 
even this relaxation technique guarantees only a relative minimum, for small de- 
viations of control point location. The sequential technique is likely to require 
more points to achieve a given level of anticipated accuracy than the minimum 
number meeting the requirements with the relaxation technique, but the sequential 
technique is, as mentioned earlier, more tractable. 

This problem, optimizing locations of control points, is analogous to the 
statistical problem of selecting a set x^ at which to evaluate the associated 
set , the combination to drive a regression between x and y. It has proven 
more practical to select an oversupply of and drop members one at a time 

according to a utility criterion, until the required performance would fail if 
any further members were dropped. A practical scheme for location of control 
points could proceed similarly. 

All too often, a desired control point location, even with the allowance of 
small deviations around it, will not contain picture structure to permit suffi- 
ciently accurate cross-correlation. The strategy for control point location 
should be robust against failure to find proper structure at a subset of the 
ideal control point locations. An iterative scheme is envisioned, in which that 
failure would send one back to the image to obtain another set of control points. 

These techniques seem clearly impractical for one-shot application. The 
amount of time spent in the computations could better have been spent in getting 
on with the registration problem by casting an overkill of control points into 
the image, and just proceeding. The application for such techniques is in a 
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production environment* where the fewest control points vith the most advantageous 
effect is desired. Using the remapping model and the a priori estimates of the 
covariance, a general optimal strategy for control point location is drawn. The 
covariance estimation is done either entirely without prejudice (E « o^I), or 
if a large number of images are to be registered against a single base image (as 
in the authors' case, multitemporal agricultural images of certain ground areas), 
that base image can be investigated for the estimations of Z. 

To the authors' knowledge, this sort of optimal placing of control points 

has not been done. The MDP uses hand-picked control points with control in any 

given frame of Landsat data coming from not only that frame but from adjacent 

4 

frames as well. In the LACIE Processor , there was essentially only a single 

point, the entire image; translational registration alone was done, with nearest- 

neighbor resampling following translational correlation done on whole-segment 

19 

basis. An experiment associated with the LACIE Processor utilized five control 

points placed at the corners and center of the image. In the DAM Package^^, manual 

control points are selected with admonishments to spread them out well in the image 

b 7 

being registered. ^ the JSC Registration Processor ' , the control points are 

placed on a uniform grid, the patch size and spacing relating so as to cause 

31 

overlap. In the "automatic” registration system for the LACIE ground truth 
segments, control points are initially placed on a uniform grid, with small 
deviations allowed in order to find a suitable correlation peak. In all the 
processors with which the authors are familiar, it has been the practice to 
follow one's instincts rather than to code in a sophisticated patch location 
algorithm. The success of those processors Indicates the eff' ” of the less- 
elaborate methods, but it remains possible to make algorithmic . i*'.rovements to 
permit sufficient accuracy with minimum computation load. 


Non-Controi Point Methods 


The previous sections dealt with the step-wise manner of first matching 
images within local patches or control points and then fitting the full mapping 
function to the resultant array of point-shifts. Although the complexity of 
straightforward multidimensional correlation (say, la translation, rotation, 
scale change, etc.) generally makes it prohibitive, a fortuitous property of 
affine distortions u^der Fourier transformation has given rise to such a method 
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which may show promise* It is shown that even with affine distortions (i.e., 
rotation^ scale change, and skew or shear distortion, in addition to pure trans- 
lation) the translation maps into a phase factor in the frequency domain, and 
the modulus of the transform function exhibits a similar, affine distortion. 

This fact, coupled with the fact that typical agricultural scenery shows con- 
siderable spatial structure, facilxtiites the use of Fourier transforms for match- 
ing in an easier way than with the images themselves. The structure of agricul- 
tural scenes generally casts the major portion of the spectral energy (i.e., 
transform modulus) along two nearly-perpendlcular axes, corresponding to the 
presence of rectangular field and road structure. In that case the axes of the 
two transformed images can be sought out by searching for lines of maxima, and 
the matching need only occur over the axes, rather than over the whole domain. 

The rotation and skew are determined by rotating corresponding axes into coin- 
cidence, and scale changes are determined by matching energy distributions along 
the axes. Once the linear portion of the affine transformation has been deter- 
mined (wholly from the modulus), the translation is determined from the phase. 

32 

The method was used successfully to determine rotation and skew distor- 
tion in airborne scanner data. It is not certain that modulus-matching along 
the axes for scale change determination is sensitive enough to give a good scale 
solution because the modulus probably changes far more slowly along an axis than 
traverse to it. This question should be addressed and the method tried on other 
data types and scene types before furthi ^ conclusions are drawn. At present, 
the method does seem applicable to agricultural images in which only transla- 
tion, rotation, and skew are present. 

Summary and Conclusions 

A number of image correlation, match offset determination, control point 
placement, and geometrical modeling techniques are in use today, but there is 
still considerable ground to cover. There appears to have been little effort 
in interchanging methods, e.g., trying different offset determination techniques 
with each matching technique, etc. Perhaps the greatest potential for correla- 
tion improvements lies in increased computational efficiency and speed. Several 
new matching techniques appear to show promise, but they need to be tested on a 
wider range of image data. There is room for improvement in the robustness of 
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offset determination schemes « The various schemes need to be compared when used 
with representative data. Improved techniques for identifying false fixes with- 
out rejecting good fixes are seriously needed. 

Methods of modeling sensor anomalies such as late line start perturbations 
need further investigation. Is the current MSS scan uonuniformity model adequate, 
or does it change with time? More attention to the degrading effects of geometric 
distortion on cross-correlation is needed by evaluating tradeoffs between patch 
size and iterative registration. Control point placement techniques need to be 
compared in terms of overall performance; particularly, the adaptive placement 
scheme described earlier needs to be tested. Robust estimation techniques saould 
be Investigated for their applicability to registration. 
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TEMPLATE MATCHING: 
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Figure A. Normalization Options For 

Registration Cross-Correlation 
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